ODYSEA Swath Simulator Examples¶
This notebook contains examples for constructing an ODYSEA orbital swath, co-locating model data, and adding expected instrument errors.
These simplified simulations are part of a much larger WaCM simulation capability at JPL, where the actual instrument, sampling, and science performance are generated.
Want to contrubite? Please do! Ask wineteer@jpl.nasa.gov how to get involved.
Imports and classes¶
In [1]:
%load_ext autoreload
%autoreload 2
In [2]:
from odysim.swath_sampling import OdyseaSwath
from odysim.errors import OdyseaErrors
from odysim.colocate_model import GriddedModel, addTimeDim
from odysim.utils import makePlot
import numpy as np
import scipy
import matplotlib
import matplotlib.pyplot as plt
import datetime
import itertools
font = {'weight' : 'bold',
'size' : 16}
matplotlib.rc('font', **font)
matplotlib.rc('font', **font)
matplotlib.rc('lines',linewidth=4)
matplotlib.rc('text',usetex=False)
matplotlib.rcParams.update({"axes.grid" : True, "grid.color": "black"})
Initialize the ODYSEA swath object¶
Orbit_fname contains an orbit propagated using the orbit height in the filename. config_fname most importantly contains the swath width to be used for swath generation. A swath is generated and centered on the propagated orbit nadir point.
In [3]:
odysea = OdyseaSwath(orbit_fname='../odysim/orbit_files/orbit_out_590km_2020_2023.npz',
config_fname='../odysim/wacm_sampling_config.py')
Create an orbit generator object between two dates¶
Note that the provided orbit is currently 2020-2023 at 590 km height in a 4AM/PM sun-synchronous orbit.
In [21]:
start_time = datetime.datetime.strptime('2020-01-22:12','%Y-%m-%d:%H')
end_time = datetime.datetime.strptime('2020-01-25:12','%Y-%m-%d:%H')
orbits = odysea.getOrbits(start_time = start_time,end_time = end_time)
Grab the first xarray dataset orbit¶
The orbit iterator we created above contains all orbits that exist between the two chosen start/end dates. Actual orbital swath data is generated on the fly whenever next() is called. A for/while loop can also be used to iterate through each orbit inside the orbits iterator. Each orbit swath will be returned as an xarray dataset containing information about sample time/lat/lon and antenna pointing information.
In [6]:
%%time
o = next(orbits)
display(o)
<xarray.Dataset>
Dimensions: (along_track: 8060, cross_track: 300)
Coordinates:
* along_track (along_track) int64 0 1 2 3 4 5 ... 8055 8056 8057 8058 8059
* cross_track (cross_track) int64 0 1 2 3 4 5 ... 294 295 296 297 298 299
Data variables:
sample_time (along_track, cross_track) datetime64[ns] 2020-01-22T20:5...
lat (along_track, cross_track) float32 79.69 79.72 ... 75.31
lon (along_track, cross_track) float32 176.9 176.8 ... 112.2
swath_blanking (cross_track) bool True True True True ... True True True
encoder_fore (along_track, cross_track) float64 -90.0 -83.37 ... 90.0
encoder_aft (along_track, cross_track) float64 -90.0 -96.63 ... 90.0
azimuth_fore (along_track, cross_track) float64 -90.0 -83.37 ... 24.55
azimuth_aft (along_track, cross_track) float64 -90.0 -96.63 ... 24.55
bearing (along_track) float64 0.0 -62.2 -77.51 ... -65.45 -65.45
Attributes: (12/14)
title: Odysea Simple Orbit Sampling V0.1
project: Odysea
summary: Simplified orbit sampling assuming basic Odysea orb...
references: Rodriguez 2018, Wineteer 2020
institution: Jet Propulsion Laboratory (JPL)
creator_name: Alexander Wineteer
... ...
geospatial_lat_min: -89.99N
geospatial_lat_max: 89.99N
geospatial_lon_min: -180.00E
geospatial_lon_max: 180.00E
time_coverage_start: 2020-01-22 20:54:00
time_coverage_end: 2020-01-22 22:29:59.640970CPU times: user 12.7 s, sys: 20.4 ms, total: 12.7 s Wall time: 12.7 s
Create a gridded model object to hold model data¶
The included GriddedModel class is designed to take in lat/lon gridded data and is used to interpolate to orbital swath data. Modify the GriddedModel class as neccessary to support your own model dataset. Much more documentation on how this works is in the GriddedModel code.
In [5]:
model = GriddedModel(preprocess=addTimeDim,n_files=200)
display(model.U)
<xarray.Dataset>
Dimensions: (time: 200, lat: 4501, lon: 9000)
Coordinates:
* time (time) datetime64[ns] 2020-01-19T21:00:00 ... 2020-01-28T04:00:00
* lat (lat) float64 -90.0 -89.96 -89.92 -89.88 ... 89.88 89.92 89.96 90.0
* lon (lon) float64 -180.0 -180.0 -179.9 -179.9 ... 179.9 179.9 180.0
Data variables:
U (time, lat, lon) float64 dask.array<chunksize=(1, 4501, 9000), meta=np.ndarray>Co-locate model currents and winds to the orbit swath¶
This code will add numerous wind and current variables into the orbit swath object. u_model, v_model, u10_model, v10_model, tx_model, ty_mode, wind_speed_model, wind_dir_model. Tx is wind stress. U10 is 10 meter wind. u is currents.
In [7]:
o = model.colocateSwathCurrents(o)
o = model.colocateSwathWinds(o)
display(o)
<xarray.Dataset>
Dimensions: (along_track: 8060, cross_track: 300)
Coordinates:
* along_track (along_track) int64 0 1 2 3 4 ... 8055 8056 8057 8058 8059
* cross_track (cross_track) int64 0 1 2 3 4 5 ... 295 296 297 298 299
Data variables: (12/17)
sample_time (along_track, cross_track) datetime64[ns] 2020-01-22T20...
lat (along_track, cross_track) float32 79.69 79.72 ... 75.31
lon (along_track, cross_track) float32 176.9 176.8 ... 112.2
swath_blanking (cross_track) bool True True True True ... True True True
encoder_fore (along_track, cross_track) float64 -90.0 -83.37 ... 90.0
encoder_aft (along_track, cross_track) float64 -90.0 -96.63 ... 90.0
... ...
u10_model (along_track, cross_track) float64 -2.985 -3.081 ... nan
v10_model (along_track, cross_track) float64 2.092 2.075 ... nan nan
tx_model (along_track, cross_track) float64 -0.01487 ... nan
ty_model (along_track, cross_track) float64 0.01042 0.01053 ... nan
wind_speed_model (along_track, cross_track) float64 3.645 3.714 ... nan nan
wind_dir_model (along_track, cross_track) float64 -54.98 -56.04 ... nan
Attributes: (12/14)
title: Odysea Simple Orbit Sampling V0.1
project: Odysea
summary: Simplified orbit sampling assuming basic Odysea orb...
references: Rodriguez 2018, Wineteer 2020
institution: Jet Propulsion Laboratory (JPL)
creator_name: Alexander Wineteer
... ...
geospatial_lat_min: -89.99N
geospatial_lat_max: 89.99N
geospatial_lon_min: -180.00E
geospatial_lon_max: 180.00E
time_coverage_start: 2020-01-22 20:54:00
time_coverage_end: 2020-01-22 22:29:59.640970Look at these model data in the swath¶
In [8]:
o.u_model.T.plot(figsize=(20,2),vmax=1,vmin=-1,cmap='RdBu_r')
o.v_model.T.plot(figsize=(20,2),vmax=1,vmin=-1,cmap='RdBu_r')
o.wind_speed_model.T.plot(figsize=(20,2))
Out[8]:
<matplotlib.collections.QuadMesh at 0x7f6ca065d6c0>
Look at these data on a geographic map¶
In [9]:
lons = []
lats = []
var_to_plot = []
# stack up some data from 2 orbits
for o in itertools.islice(orbits, 2):
o = model.colocateSwathCurrents(o)
lons = np.append(lons,o.lon)
lats = np.append(lats,o.lat)
var_to_plot = np.append(var_to_plot,o.u_model)
# make a map with utils.makePlot(); this will grid data for us.
makePlot(lons,lats,var_to_plot,vmin=-1,vmax=1,cblabel='U Current [m/s]',colormap='RdBu_r',figsize=(20,10),bg=True,gridMe=True)
Initialize an OdyseaErrors class that will do uncertainty estimation.¶
In [11]:
e = OdyseaErrors(type='cbe') # 'cbe' (current best estimate) or 'mev' (maximum expected value, which includes 3dB SNR margin)
Set surface current uncertainties (radial and geographic) using a simple constant wind speed error model.¶
This will set all expected standard deviations and resulting uncertainties. Using simulated_baseline will generate uncertainties that are a function of cross-track location and wind speed. If winds have been co-located to the swath from a model, those wind speeds and directions will be used to generate the error field, otherwise, a constant wind speed/dir will be used from the function arguments.
In [12]:
o = next(orbits)
o = e.setCurrentErrors(o,etype='simulated_baseline',wind_speed=7.5,wind_dir=0)
display(o)
Using wind speed from function args. Assign "wind_speed" and "wind_dir" variables to orbit dataset if desired.
<xarray.Dataset>
Dimensions: (along_track: 7976, cross_track: 300)
Coordinates:
* along_track (along_track) int64 0 1 2 3 4 5 ... 7971 7972 7973 7974 7975
* cross_track (cross_track) int64 0 1 2 3 4 5 ... 294 295 296 297 298 299
Data variables: (12/17)
sample_time (along_track, cross_track) datetime64[ns] 2020-01-23T01:4...
lat (along_track, cross_track) float32 79.79 79.82 ... 75.08
lon (along_track, cross_track) float32 83.85 83.66 ... 43.16
swath_blanking (cross_track) bool True True True True ... True True True
encoder_fore (along_track, cross_track) float64 -90.0 -83.37 ... 90.0
encoder_aft (along_track, cross_track) float64 -90.0 -96.63 ... 90.0
... ...
vx_std (along_track, cross_track) float64 1.145e+15 ... 2.059e+15
vy_std (along_track, cross_track) float64 -0.07014 ... 0.1261
u_std (along_track, cross_track) float64 8.403e+07 ... 1.786e+15
v_std (along_track, cross_track) float64 1.145e+15 ... 1.025e+15
u_error (along_track, cross_track) float64 -4.157e+07 ... 7.477e+14
v_error (along_track, cross_track) float64 2.857e+15 ... 2.713e+14
Attributes: (12/14)
title: Odysea Simple Orbit Sampling V0.1
project: Odysea
summary: Simplified orbit sampling assuming basic Odysea orb...
references: Rodriguez 2018, Wineteer 2020
institution: Jet Propulsion Laboratory (JPL)
creator_name: Alexander Wineteer
... ...
geospatial_lat_min: -89.99N
geospatial_lat_max: 89.99N
geospatial_lon_min: -180.00E
geospatial_lon_max: 180.00E
time_coverage_start: 2020-01-23 01:44:00
time_coverage_end: 2020-01-23 03:18:59.542614Show the geographic u/v surface current uncertainties.¶
In [13]:
o.u_error.T.plot(figsize=(20,2),vmax=.5,vmin=-.5,cmap='RdBu_r')
o.v_error.T.plot(figsize=(20,2),vmax=.5,vmin=-.5,cmap='RdBu_r')
Out[13]:
<matplotlib.collections.QuadMesh at 0x7f6cc0297af0>
Repeat the error estimation, but this time, co-locate model winds to the swath and use them to set uncertainties.¶
In [14]:
o = model.colocateSwathCurrents(o)
o = model.colocateSwathWinds(o)
o = e.setCurrentErrors(o,etype='simulated_baseline')
display(o)
<xarray.Dataset>
Dimensions: (along_track: 7976, cross_track: 300)
Coordinates:
* along_track (along_track) int64 0 1 2 3 4 ... 7971 7972 7973 7974 7975
* cross_track (cross_track) int64 0 1 2 3 4 5 ... 295 296 297 298 299
Data variables: (12/25)
sample_time (along_track, cross_track) datetime64[ns] 2020-01-23T01...
lat (along_track, cross_track) float32 79.79 79.82 ... 75.08
lon (along_track, cross_track) float32 83.85 83.66 ... 43.16
swath_blanking (cross_track) bool True True True True ... True True True
encoder_fore (along_track, cross_track) float64 -90.0 -83.37 ... 90.0
encoder_aft (along_track, cross_track) float64 -90.0 -96.63 ... 90.0
... ...
u10_model (along_track, cross_track) float64 -2.764 ... -5.519
v10_model (along_track, cross_track) float64 0.09098 ... -7.488
tx_model (along_track, cross_track) float64 -0.01045 ... -0.07015
ty_model (along_track, cross_track) float64 0.0003438 ... -0.09517
wind_speed_model (along_track, cross_track) float64 2.766 2.884 ... 9.302
wind_dir_model (along_track, cross_track) float64 -88.11 ... -143.6
Attributes: (12/14)
title: Odysea Simple Orbit Sampling V0.1
project: Odysea
summary: Simplified orbit sampling assuming basic Odysea orb...
references: Rodriguez 2018, Wineteer 2020
institution: Jet Propulsion Laboratory (JPL)
creator_name: Alexander Wineteer
... ...
geospatial_lat_min: -89.99N
geospatial_lat_max: 89.99N
geospatial_lon_min: -180.00E
geospatial_lon_max: 180.00E
time_coverage_start: 2020-01-23 01:44:00
time_coverage_end: 2020-01-23 03:18:59.542614Show the geographic u/v surface current uncertainties.¶
In [15]:
o.u_error.T.plot(figsize=(20,2),vmax=.5)
o.v_error.T.plot(figsize=(20,2),vmax=.5)
Out[15]:
<matplotlib.collections.QuadMesh at 0x7f6b407127a0>
In [18]:
lons = []
lats = []
var_to_plot = []
# stack up some noisy current data from 2 orbits
for o in itertools.islice(orbits, 2):
o = model.colocateSwathCurrents(o)
o = model.colocateSwathWinds(o)
o = e.setCurrentErrors(o,etype='simulated_baseline')
lons = np.append(lons,o.lon)
lats = np.append(lats,o.lat)
var_to_plot = np.append(var_to_plot,o.u_model+o.u_error)
# make a map with utils.makePlot(); this will grid data for us.
makePlot(lons,lats,var_to_plot,vmin=-1,vmax=1,cblabel='U Current [m/s]',colormap='RdBu_r',figsize=(20,10),bg=True,gridMe=True)
Plot ~two days worth of orbits. Do not attempt to remove any regions of large error.¶
In [22]:
lons = []
lats = []
var_to_plot = []
# stack up some noisy current data from 2 orbits
for o in itertools.islice(orbits, 30):
o = model.colocateSwathCurrents(o)
o = model.colocateSwathWinds(o)
o = e.setCurrentErrors(o,etype='simulated_baseline')
mask = np.abs(o.u_error) < 9999
lons = np.append(lons,o.lon[mask])
lats = np.append(lats,o.lat[mask])
var_to_plot = np.append(var_to_plot,o.u_model[mask]+o.u_error[mask])
# make a map with utils.makePlot(); this will grid data for us.
makePlot(lons,lats,var_to_plot,vmin=-1,vmax=1,cblabel='U Current [m/s]',colormap='RdBu_r',figsize=(20,10),bg=True,gridMe=True)
Mask out non-physically large currents (large due to error). Note this could be done by looking at the error itself earlier.¶
In [25]:
mask = np.isfinite(var_to_plot) & (np.abs(var_to_plot)<2) # assumes total current must be less than 2 m/s to mask out very large errors
makePlot(lons[mask],lats[mask],var_to_plot[mask],vmin=-1,vmax=1,cblabel='U Current [m/s]',colormap='RdBu_r',figsize=(20,10),bg=True,gridMe=True)
In [ ]: